Anharmonicity in one-dimensional electron-phonon system 
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We investigate the effect of anharmonicity on the one-dimensional half-filled Holstein model by 
using the determinant quantum Monte Carlo method. By calculating the order parameters we 
find that with and without anharmonicity there is always an transition from a disorder phase to 
a dimerized phase. Moreover, in the dimerized phase a lattice dimerization and a charge density 
wave coexist. The anharmonicity represented by the quartic term suppresses the dimerization as 
well as the charge density wave, while a double-well potential favors the dimerization. In addition, 
by calculating the correlation exponents we show that the disorder phase is metallic with gapless 
charge excitations and gapful spin excitations while in the dimerized phase both excitations are 
gapful. 



I. INTRODUCTION 

The interaction between electrons and ion vibra- 
tions is responsible for many fundamental phenomena 
in solids. For example, as shown in BCS theory in 
electron-phonon(EP) interaction induces a weak attrac- 
tion between electrons near Fermi surface with oppo- 
site momenta and opposite spins and eventually leads to 
superconductivity!. In quasi one-dimensional materials 
it is well known that EP coupling usually causes metal- 
insulator transition. In the insulating phase a band gap 
opens at the Fermi energy and simultaneously lattice is 
distorted resulting in a larger unit cell. This transition 
is called Peierls transition^. 

For simplicity most of theoretical works consider only 
harmonic ion vibration and in only a few works higher 
orders are taken into account. For example, the effect 
of anharmonicity on superconductivity was studied by 
Freericks et al.'^ more than ten years ago in infinite- 
dimensional limit. However, it seems that the anhar- 
monicity represented by the quartic term does not en- 
hance the superconductivity transition temperature as 
expected. Another example is that polaronic properties 
with anharmonicity were investigated by Chatterjee and 
Takada'* some time ago. 

Recently a renewed interest in anharmonic ion vi- 
bration was raised in /3-pyrochlore oxides, especially in 
KOS2O6. Experiments show that it has unusual behav- 
iors in both normal state and superconductivity stated""—. 
In KOS2O6 an Os-0 network forms a oversized cage 
and K ion sits in the cage. As the mass of K ion is 
small compared with Rb or Cs it oscillates relatively far 
from the equilibrium position. Density functional the- 
ory shows that this movement is highly anharmonic and 
the harmonic term of the potential can even be zero 
or negativeiSiiii. Some of the anomalous behaviors ob- 
served in KOS2O6 was explained by Dahm and one of the 
present authors by treating the anharmonic ion vibration 
in a quasi-harmonic approximation with temperature- 
dependent frequencyji^ii^. The spectral function and 
NMR relaxation rate have been discussed for a model 



with a general anharmonic potential^*. Otsuka et al^^ 
has consider the effect of static anharmonic potential 
on quasi-one-dimensional extended Hubbard model and 
rich phase diagrams have been found at finite tem- 
perature. However, in quasi-one-dimensional materials, 
phonon fluctuation is usually important^^. 

In this paper we will use the one-dimensional spinfuU 
Holstein model as a prototypical model to study the ef- 
fects of anharmonicity of ion vibrations. Since here we 
are more interested in the anharmonicity the model we 
will study is given by 
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(1) 



where pi and are the momentum and position operator 

of ions at site i, V(rir2 ■ • • r^) = J2i=iiT''"i ~^ '^'''t^ 
potential energy of ion vibration with L being the lattice 
size. 7 is always nonnegative so that the ion potential is 
bounded, t is the hopping term of electrons onto nearest- 
neighbor sites. (T is the spin index and in our work we 
only consider spinful model and then a takes two values 
a =tii- cj^lcicr) is the creation(annihilation) operator 
for an electron at site i with spin index a, riicr — c\^Cia 
is the electron number operator at site i for spin tr. The 
electrons couple locally with ions and the coupling con- 
stant is given by A. /i is the chemical potential for con- 
trolling the electron number. If = Hamiltonian 
([T]) is invariant under the particle-hole transformation 
Pi —Pi,ri — ?> —ri,Ci — ?• (— l)'cj and thus it is always 
half- filled at pi — 0. 

If 7 = Hamiltonian ([I} reduces to the standard Hol- 
stein model. In the past several decades it has been exten- 
sively studied by many authors. The mean-field theory 
and the expansion from the strong coupling limit predict 
that at half-filling for spinless case there is a phase transi- 
tion from a disorder phase to a dimerized phase while for 
spinful case the ground state is dimerized for any finite 
EP coupling and any finite phonon frequencj»i^. How- 
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ever, the latter conclusion was challenged by Wu et aU^. 
Based on functional integral analysis they argued that 
the fluctuation of phonons could also destroy the dimer- 
ization for spinful case when A is smaller than a critial Ac. 
This conclusion was confirmed later by accurate density- 
matrix renormalization group calculations^ as well as 
quantum Monte carlo simulations^^. 

However, it is not clear yet what is the effect of anhar- 
monicity 7 7^ in one dimension. In particular, if the 
harmonic term K is zero or even negative, the ion po- 
tential becomes double-well. Then a natural question 
is whether there are some fundamental changes com- 
pared with the harmonic case? In the following we 
will clarify these issues by investegating the Hamilto- 
nian ([T]) by using the determinant quantum Monte Carlo 
simulationa^^i^^ . In this model the electrons of two spin 
species couple equally to ions, so there is no sign problem. 
In our simulations the largest lattice size is up to 32 and 
finite-size scaling is used to extrapolate our results to the 
thermodynamic limit. The inverse temperature j3 is up 
to 24. We check our results with different /3 to confirm 
that the temperature in our simulation is low enough to 
extract ground state properties. The time slice Ar = 0.1 
is used in most of our simulations and we also check the 
results carefully by calculations with smaller At and con- 
firm that the results are not sensitive to At. 

In this paper we consider only the half-filled case by 
putting /i = 0. Moreover for simplicity we also fix 
m = 0.5 and < = 1. In Section II, firstly we will take 
io — ^ Kim = 1, and show the results as a function 
of A. We find that for both 7 = and finite 7 there 
is a phase transition from the disordered phase to the 
dimerized phase. The dimerization occurs accompanied 
by the charge density wave ( CD W) transition simultane- 
ously. Secondly we show the phase diagram in (A, 7) 
plane. Finally we show the order parameters nip and 
TOe as a function of K with fixed A = 1.4 and 7 — 0.1. 
Similarly, a transition from the dimerized CDW phase 
to the disordered phase is found. In Section III, we will 
present results on static charge structure factor and spin 
structure factor. These data confirm our previous conclu- 
sions. By calculating correlation exponents, we find that 
for any finite A there is a spin gap. However, a charge 
gap opens only after A > Ac, or < K^. In Section IV, 
we summarize our results. 
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FIG. 1: (Color online)Finite-size scaling of square of the order 
parameters mp(open symbols) and mg (filled symbols). The 
parameters (A, 7) for black circles, red squares, green dia- 
monds and blue triangles are (0.7, 0.0), (1.0, 0.1), (1.0, 0.0) and 
(2.0,0.1), respectively. 

change as a function of A when the anharmonic term 7 
is introduced. For our purpose we define the staggered 
correlation function of lattice displacements 
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and the staggered correlation function of the electron 
densities 



(3) 



The order parameters rrip and me are correspondingly 
defined as 



and 
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II. DIMERIZATION AND CHARGE DENSITY 
WAVE ORDER PARAMETERS 

Since Hamiltonian ([T]) with 7 = has been extensively 
studied then we will focus on the case with a finite 7 in 
this work. Data for 7 = are calculated just for compar- 
ison with the DMRG results^S. When 7 is larger than 
zero the anharmonic term increases the elastic energy 
of ion. Therefore a straightforward conjecture is that it 
will suppress the dimerization as well as the CDW order. 
In the first step we examine how the order parameters 



In the thermodynamic limit is expected to be finite in 
the dimerized phase while it is zero in the disorder phase. 
Similarly, TOc is finite in the CDW phase and zero in the 
disorder phase. To be specific we first fix = 0.5. In 
Fig. [l]we show both and as a function of inverse 
length of chain with various (A, 7). Now let us see the 
data without anharmonicity, i.e., 7 = 0. As shown by 
circles in (a) both and clearly extrapolate to zero 
in the thermodynamic limit for A — 0.7 while for A = 
1.0(shown by diamonds in (b)) the order parameters are 
finite. This fact tells clearly that there is at leaf one phase 
transition at Ac satisfying 0.7 < Ac < 1.0. Next we turn 
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FIG. 2: mp(open symbols) and (filled symbols) as a func- 
tion of A for 7 = O(left) and 7 = 0.1 (right). The crosses in 
the left panel are me x X/K. 

to see the data with 7 > 0. As a representative example, 
here we choose 7 = 0.1. As A = 1.0(shown by squares 
in (a)) both and vanish in the thermodynamic 
limit, which is different from 7 = case. However, as 
we increase EP coupling A we reenter in the dimerized 
CDW phase, which is shown by triangles in Fig. [IJb) 
with A = 2. As we have seen, a larger Ac is required when 
the anharmonicity 7 is present and this is consistent with 
our conjecture. 

To extract the critical point Ac accurately we calcu- 
late the order parameters by scaning A with fixed 7. The 
order parameters as a function of A are shown for both 
7 = in the left panel of Fig. [2]and 7 = 0.1 in right panel 
of Fig. [51 respectively. For 7 = the data obtained by 
the QMC are in good agreement with those obtained by 
the DMRG^". The common feature for both 7 = and 
7 — 0.1 is that both rrip and me vanish at the same 
critical points, and thus there is one critical point. The 
critical values Ac are estimated to be 0.8 and 1.4, respec- 
tively. For 7 = 0, rup and me is found numerically^Sito 
satisfy 
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(6) 



However, this relation does not hold for 7 = 0.1. In 
particular, rrip/me is smaller than X/K. It is plausible 
to say that the anharmonic potential 7 has "stronger" 
suppression on the ion dimerization than on the CDW 
formation. In Fig. [3] we show the schematic phase di- 
agram in (A, 7) plane for K = 0.5. We find that as 7 
increases the critical value A also increases. This results 
also support our conclusion that the anharmonicity sup- 
presses the dimerization and CDW order. 

As the coefficient of the harmonic term of the poten- 
tial may be negative and then the total ion potential has 
a form of a double-well, we also calculate the order pa- 
rameters as a function of K with A = 1.4 and 7 = 0.1, 
which are shown in Fig. In this figure we see that 
when K < Kc = 0.5 the ground state is in the dimerized 
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FIG. 3: (Color online) Phase diagram in (A, 7) plane for K = 
0.5. The solid line is a guide of eyes. 

1.5 



A~Am„ 



1=IA, y=0.\ 




FIG. 4: (Color online)Order parameters rUp and rrie are shown 
as a function of K. Here A is set to be 1.4 and 7 = 0.1. 

CDW phase while for K > Kc it is in the disorder phase. 
No new phases are found in our simulations. Thus we 
may conclude that negative K favors the dimerization 
with the CDW order but does not introduce new phases 
in this model in one dimension. 



III. STATIC STRUCTURE FACTOR AND 
CORRELATION EXPONENTS 

One-dimensional correlated electronic systems 
can usually be described by Tomonaga-Luttinger 
liquids (TLL)^. In the TLL theory charge degree of 
freedoms and spin degree of freedoms are seperated. 
In particular, decay of the charge and spin correlation 
functions are determined by two exponents Kp and K„^ 
respectively. To explore properties of the disorder phase 
and the dimerized one further we calculate the charge 
structure factor and spin structure factor, which are 
defined by 
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FIG. 5: (Color online) Top: Static charge structure factor 
C{q) for A = 1.0. Inset: C{%)/L is shown as a function of l/i. 
Bottom: Static charge structure factor for A = 2.0. Inset: 
C{ti)/L is shown as a function of 1/L. 



and 



In particular, Kp and K^, can be calculate d^ by 
Kp = TT lim dC{q)/dq. 



and 



Ka — 47r lim dS{q)/dq, 



(8) 



(9) 



(10) 



In Fig. [5] we show C{q) as a function of q for A = 1.0 
and A — 2.0 with K = 0.5 and 7 — 0.1. In the top panel 
we find that C{q) for L — 16,24 and 32 fall almost into 
one curve and a peak is present at g = tt. In the in- 
set we show C{'k)/L as a function of 1/L. It vanishes in 
the thermodynamic limit. This figure demonstrates that 
our sizes are large enough and finite size effect is negligi- 
ble. One can conclude that there is only a dominant 2kF 
short-range correlation and no long-range order. In the 
bottom panel it is interesting to notice that C{q) exhibits 
a singularity at 5 = tt and C{tt) /L is finite in the thermo- 
dynamic limit. This figure indicates that a long-range 
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FIG. 6: Static spin structure factors for A = 1.0 and A = 2.0 
with the lattice size L — 16, 24 and 32. K is chosen to be 
0.5 and 7 is 0.1 and thus the parameter set used in the top 
panel is in the disorder phase and that of the bottom is in the 
dimerized phase. 

charge-charge correlation is present. These two figures 
are consistent with our conclusion that A = 1 is in the 
disordered phase and A = 2 is in the dimerized CDW 
phase. 

Now we turn to study of the spin correlation functions. 
In Fig. |6]we show the data of static spin structure factor 
for various lattice sizes for A — 1.0 and 2.0 at K = 0.5 
and 7 = 0.1. First, from the collapse of numerical data 
of various lattice sizes we may conclude that finite size 
effect is negalible in the data. Secondly we observe that 
a peak is present at q — n for both figures. Moreover, the 
height of the peak is found to be finite in both A = 1.0 
and A = 2.0. We then conclude that in this model there 
is a dominant short-range 2kF spin correlation in both 
the disorder phase and the dimerized phase. 

Next we will analyse the charge and spin correlation 
exponents Kp and K^. Since our model is invariant un- 
der spin rotation, we expect that if it is in the TLL phase 
is equal to 1. However, if it is less than 1 a spin gap 
is expected to oper*2iiS^. This quantity is quite sensitive 
to the opening of a spin gap and therefore we use it as 
a criterion to determine whether there is a spin gap or 
not. In the top panel of Fig. [7] we show both Kp and 
Ka as a function of A at 7 = 0.1. The behaviors are ac- 
tually qualitatively similar to the results of the standard 
Holstein model, which are obtained by the Monte carlo 
simulationsSi. Namely, in the weak coupling region, Kp 
is larger than 1. Moreover, it shows nonmonotonous de- 
pendence on A: it increases from 1 as A increases from 
zero, and after it reaches a maximum it starts to decrease. 
In the strong coupling region Kp is always smaller than 
one. The crossover point with the line Kp = 1 is esti- 
mated to be between A = 1.3 and A = 1.4. This value is 
in agreement with the critical point Ac — 1.4. We believe 
the small deviation is due to numerical errors. In the 
weak coupling region Kp > 1 is intepreted as the attrac- 
tion of effective electron-electron interaction mediated by 
phonons while in the strong coupling region Kp < 1 is 
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FIG. 7: (Color online) Top: Kp and Ka as a function of A for 
K = 0.5 and 7 = 0.1. Bottom: Kp and Ka as a function of 
is: for A = 1.4 and 7 = 0.1. 

intepreted as effective repulsive electron-electron inter- 
action in the dimerized phase. However, Kcr does not 
show such crossover and it always smaller than 1 and 
monotonously decreases as a function of A. We then con- 
clude that a spin gap exists for any finite A. In the bottom 
panel we show Kp and K^ as a function of -RT for 7 = 0.1. 
In this figure negative K is also considered. Again Kp 
crosses from below 1 to above 1 as if increases. The 
cross point seemingly is consistent with the critical point 
Kc = 0.5 within our error bars. In the same way as be- 
fore, Ka is small than 1 in the whole range, indicating the 
presence of a spin gap. Compared with the results ob- 



tained by the DMRG^° and the QMC^^ we may conclude 
that the disorder phase is Luther-Emery liquid. 



IV. SUMMARIES 

In summary, in this paper we have investigated the ef- 
fect of anharmonicity on the one-dimensional half-filled 
Holstein model. By studying order parameters, static 
charge and spin struture factors and correlation expo- 
nents we find a transition from Luther-Emergy liquid to 
dimerized CDW phase. This is essentially similar to stan- 
dard Holstein model. The effect of anharmonicity is to 
suppress the dimerization as well as the charge density 
wave. Recently, the effective interaction of electrons me- 
diated by anharmonic phonons has been derived and the 
transition temperature of CDW of Hamiltonian ([T]) has 
been discussed^^. The transition temperature decreases 
monotonously as 7 increases. This behavior is also con- 
sistent with the present conclusions. On the other hand, 
a double-well potential with a negative quartic term fa- 
vors the dimerization and the CDW order. 
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